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Abstract 

It has been asserted previously by the author that optimal control theory can be a valuable framework for theoretical studies 
about the shape that a wind instrument should have in order to satisfy some optimization criterion, inside a fairly general class. 
The purpose of the present work is to develop this new approach with a look at a specific criterion to be optimized. In this 
setting, the Webster horn equation is regarded as a controlled dynamical equation in the space variable. Pressure is the state, 
the control being made of two parts : one variable part, the inside diameter of the duct and one constant part, the weights 
of the elementary time-harmonic components of the velocity potential. Then one looks for a control that optimizes a criterion 
related to the definition of an oscillation regime as the cooperation of several natural modes of vibration with the excitation, 
the playing frequency being the one that maximizes the total generation of energy, as exposed by A.H. Benade, following H. 
Bouasse. At the same time the relevance of this criterion is questionned with the simulation results. 

^ _ I. Introduction 

ps) ' Designing high quality musical instruments has been the main concern of makers for centuries, accumulating know-how 
from their predecessors by some kind of trial-and-error process, involving musicians. On another hand, since Bernouilli and 
Lagrange in the XVIII*'' century and the important writings of Helmholtz in the XIX*'* century a great deal of research 
'"s has also been conducted in order to understand the principles underlying sound production in traditional or modern musical 
On instruments : see e.g. [2], [3], [6], [7], [13] for presentations of this vast subject. Wind instruments such as brasses and 
woodwinds, making use of an air column, must be designed in such a way as to properly arrange the natural frequencies 
of this air column, in order that regimes of oscillation can take place in conjunction with the nonlinear flow-control device 
(reed or player's lips) [2]. As a consequence, harmonicity requirements for wind instruments between natural frequencies 
are desirable properties, due to the fact that intonation, responsiveness and tonal colour can be attributed to well established 
physical properties of the instrument [1]. It is generally accepted that these harmonicity requirements can lead to rather 
1^ different shapes, which can be (piecewise) cones, cylinders or more complex, as can can be seen from the actual duct of 
real wind instruments. Nevertheless, see [9] for interesting results about piecewise conic ducts, following an approach based 
^ on characteristic impedances and transfer matrices. Thus there is room for developping new methods that could help for the 
I— I design. The present work is dedicated to such a task, with main accent on the methodological aspects. 

In [18], transmission line modelling of horns is used together with finite dimensional optimization techniques and a practical 
^ tool [19] was developed for instrument makers, in order to help improving existing instruments as well as to design new ones 
according to a given specification. Desired properties, specified in musical terms like intonation, response and pitch variability 
T— I were compared with calculated values based on an instrument's actual geometry. Similarly, in [22], a method is presented for 
optimizing the shape of a brass instrument with respect to its intonation and impedance peak magnitudes. The instrument is 
modelled using a one-dimensional transmission line analogy with truncated cones. Through the use of an appropriate choice 
T— I of design variables, the finite dimensional optimization finds smooth horn profiles, that can also help in correcting existing 
instruments, the shape having been designed from an a priori choice of a succession of cones. In [14], a frequency-domain 
method, using inverse quantum scattering for the one-dimensional Klein-Gordon equation, allows to recover the area function 
. . , of a given acoustical duct in a noninvasive way, without measuring directly neither the input impedance nor the reflectance : 
^ ' this last one is mathematically derived from the wave radiated in response to a high-impedance source. In [16], parameter 
k> optimization techniques have been used to design shapes for brass-trapping Helmholtz resonators that resonate at a design 
^ set of acoustic eigenvalues, taking into account physical and geometrical constraints. In [15], a finite-element eigenanalysis 
5^ , model of bars is used to compute optimal shapes for mallet percussion instruments such as vibraphone or marimba-type 
bars. The objective function for the optimization procedure is a targeted set of modal frequencies. 

As one can see, optimization is at the heart of research on the design question for musical instruments and the eigenstructure 
often plays a central role in the optimization criteria. The shape optimization of wind instruments, an inverse problem, appears 
to be a complex one, when compared to the analysis of what physically occurs inside a given instrument, a direct, somewhat 
easier albeit not at all obvious, problem, as can be seen from the great amount of research that it is the subject of. The 
purpose of the present work is a continuation of preliminary ones by the author [20], [21], where the design question 
for axisymmetric wind instruments was posed as an optimal control problem for the so-called Webster equation through a 
suitable reformulation as a dynamical equation in the spatial dimension, thereby using infinite dimensional optimization of 
optimal control theory for continuous systems. This is in sharp contrast with the above mentionned references which use 
finite-dimensional optimization schemes, i.e. with a finite number of parameters as unknowns in the theoretical formulation. 
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For more details about the horn equation, see [12] or the more recent [25] for a more mathematical and updated exposition. 
Here, besides the exposition of how optimal control theory fits the wind instruments design question, a special attention is 
given to the choice of one specific criterion, on the basis of physical considerations, whereas in [20], only a fairly general 
class of criteria -affine in the control- were considered that led to a general result about the structure of the shape for the 
duct of a wind instrument. 

The chosen criterion is in no way unique : many other choices could be made, depending on which characteristics (physical, 
perceptive and so on) are focused on. It appeared to be interesting to investigate from a musical acoustics point of view 
as it comes from one definition of an oscillation regime inside wind instruments, understood as the cooperation of several 
natural modes of vibration, the 'playing frequency' being the one that maximizes the total generation of energy [3], [2]. At 
this stage, it has some arbitrariness in it : the physical relevance of this criterion is not taken here for granted, as no other 
published material has been found (than that of A.H. Benade) using this notion of oscillation regime. The reason for the 
choice of this criterion is twofold : first, to illustrate the general methodology that uses optimal control theory ; second, to 
question the relevance of this criterion, that appeared to be appealing at first sight. The purpose then is to question how 
interesting it is for designing wind instruments more than to assert such a relevance a priori : numerical results in section |V] 
show some qualitative features linked to the chosen criterion. 

Nevertheless, on the basis of this criterion, in the context of optimal control theory, necessary conditions coming from 
the strong Pontryagin Maximum Principle [23], [8] are derived and some numerical simulations are conducted, in order 
to illustrate the theoretical predictions. In the presented approach, the Webster horn equation is regarded as a controlled 
dynamical system in one spatial variable (the axis of an axisymmetric instrument). The pressure (or the velocity potential) 
is regarded as the state and the control is made of two parts : one is the variable diameter (actually its derivative with 
respect to the axial space dimension) of the duct. The second part of the control is made of the constant weights of a 
time-harmonic decomposition of the velocity potential, which are unknown, according to the chosen criterion, the energy 
inside the duct. Then one looks for a control that maximizes the latter. This is another difference with the approach in [18], 
[22], where the measure is taken to be the deviation from a measured input impedance peaks location, making use of a 
model-fitting oriented approach. Whereas these works can be described as inverse problems, the present one falls into the 
category of design problems. But this distinction is a matter of convenience as both problems only differ from each other 
by the chosen criterion, the unknown being the duct shape in each case. It is worth noting in that respect that a similar 
criterion as in [22] could be chosen here as well, while defining the corresponding output of the model, leading to the same 
problem of fitting a model to given data. This will be the subject of future research, with interesting comparisons to be 
made, the main point here being methodological. The approach followed hereafter shows that, with the chosen optimization 
criterion and the modelling approach, the obtained shape is piecewise continuous. Moreover, when taking into account more 
and more time-harmonic components, is appears from the simulations that the the obtained duct shape is a cone : such a 
result tends to confirm that this criterion is likely not to be sufficient to design high-quality instruments but the realistic 
mode shapes that are obtained are encouraging. It should be noted that, to the best knowledge of the author, this way of 
looking at the classical Webster horn equation, as a controlled dynamical equation, appears to be new. One important point 
is that extensions of the optimal control methodology to other more realistic, nonlinear models of wind instruments and 
design criteria can be considered, although Fourier analysis will generally be no more valid in such contexts : the choice of 
an adequate criterion will remain crucial. Numerical simulations show that the proposed approach is appealing for design 
purposes. Nevertheless, the presented work has to be considered as a preliminary study that opens new perspectives for 
the design of wind instruments : it is not to be understood as an achievement in itself and has surely to be investigated 
more deeply for real applications. In particular, several aspects are not touched here such as taking account of tone holes in 
woodwinds (flutes, oboes, etc.), incorporating the nonlinear excitation mechanism (reed, lips for brasses) or choosing other 
criteria, in order to include perceptive parameters (such as given by a musician), etc.. .Much remains to be done. 
The paper is organized as follows : section |lll in order to fix the notations, recalls the standard horn equation and the 
definition of an oscillation regime as used here. Basic results from Optimal Control theory are briefly given and lead to 
the formulation of the horn equation as a controlled dynamical system. Then a detailed exposition of how the oscillation 
regime is formulated is given in section |lll] Section |IV] gives the mathematical formulation of the optimization problem to 
solve, through the necessary conditions coming from optimal control theory. Numerical results are given in section |V] : they 
have to be considered as preliminary results as no thorough or exhaustive investigation have been conducted at this stage. 
Eventually, some conclusions and perspectives are drawn in section |VT] An appendix provides a variational derivation of 
the horn wave equation and reminds basic results from Sturm-Liouville theory. 

II. Preliminaries : Horn wave equation. Oscillation regime. Optimal control theory 

A. Wave horn equation, oscillation regime 

Consider an axisymmetric horn with length L and section diameter D{x), a function of the independent variable, the 
space dimension x, according to Fig. 1. Let po be the medium mass density, p{x, t) the acoustic pressure in the medium, i.e. 
the deviation from the atmospheric pressure, and v{x, t) the particular velocity inside the horn at abcissa x. Let fc = w/c be 




the wave number, c being the sound velocity and u the angular frequency in case of a harmonic regime that is focused on 
here. For derivatives of any quantity the index notation will be used in this section : for example the time (resp. space) first 
derivative of cf) will be denoted 0f (resp. (f>x)- Let us first recall some usual modelling hypothesis : i) the diameter D{x) is a 
slowly varying function of x ; ii) acoustical quantities are functions of x and t only ; iii) a velocity potential cf) exists, which 
means that v — Vcf) = (px ; iv) one focuses on standing wave, time-harmonic solutions that write : 4>{x,t) = 5R((y5(a;)e*'^*) 
but one will omit the 5R part in the sequel as it is related to the time dependence that will not be involved in the optimization 
process below ; v) in this work there is no dissipation. It can be shown that, in this situation, the pressure at x writes : 
p = —po(pt- It is worth noticing that the first hypothesis is necessary for the ID plane wave, Webster equation to be valid [16] 
or [7], chap. 6. Another limitation to this equation is to lower frequency modes. In [7], the following approximate condition 
on the radius derivative R is given for the horn equation to be valid for plane waves : 

1 /2 

- kR {x)dx < 1 (1) 
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This condition will be used to fix bounds on the control in the simulations presented in section |V] With these conventions, 
the Webster horn equation is easily obtained as the Euler equation of a variational problem : 

-^4itt - 4>xx - 2—^(f>x = (2) 

which reduces to : 

<p" +2^(1)' +e (1)^0 (3) 

in case of a harmonic regime of interest here (see appendix a for details). Boundary conditions are then given for this 
equation, depending on the physical situation under study, making it a Sturm-Liouville problem (see appendix b for a brief 
reminder of two useful results). 

On another side, the notion of oscillation regime can be defined, following H. Bouasse [3] and A.H. Benade [2]. It is "that 
state of the collective motion of an air column in which a nonlinear excitation mechanism collaborates with a set of air 
column modes to maintain a steady oscillation containing several harmonically related frequency components, each with its 
own definite amplitude. Then, the 'playing frequency' is the one that maximizes the total generation of energy" [2]. 
Chosing such a definition as the basis of an optimization criterion relies upon physical/acoustical considerations of what 
happens inside the instrument. Nevertheless, it has some arbitrariness when considering perceptual quality of sound, as 
discussed in the introduction. But further investigations looked desirable. Other criteria will be considered in future work 
so the present choice is only a first step towards physical/perceptual accounting in wind instruments design. 
Notice also that the complex question of the nonlinear excitation mechanism, as the reed of an oboe or the lips of a brass 
player e.g., is discarded (equation (O is homogeneous), although nonlinear effects can be important due to the nature of 
fluid dynamics, e.g. inside a reed or at the excitation localization in a flute. The hypothesis is made here that the nonlinear 
effects are, in a first approximation, encoded in the Fourier coefficients of the potential (see section HIH). corresponding to an 
arbitrary and unknown excitation. This is surely not satisfactory from a physical point of view and will have to be modified 
carefully in the future. 

In the present context, where losses are assumed to vanish or to be included into the boundary conditions, the energy 
inside the duct, noted E, is a conserved quantity : 

E^nx,t) + U{x,t)^^-f f D^{\4>,\^ + \'k\^)dx (4) 
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with ^ — and T{x, t), U {x, t) are computed in the appendix. Following the above definition of an oscillation regime, the 
energy can be taken as a design criterion for the duct shape. As it is a function of the diameter, this one can be considered as 



an unknown function, provided convenient data are given for the design. In that respect, as musical instruments are designed 
to operate in harmonic regimes, the potential, being in that case periodic, is amenable to a Fourier series analysis and can 
thus be written as a linear combination of time-harmonic elementary components, the coefficients of which will be other 
unknowns of the design problems, while the frequencies will be considered as the data. This will be made more precise in 
section |lll] Thanks to the computations in the appendix, the link is established between an oscillation regime in the sense of 
A.H. Benade and the way the Webster horn is derived, as both rely upon the same definite integral, namely the lagrangian 
action, through a Legendre transformation. As optimal control theory is variational too, one has a coherent set of tools for 
the design methodology to be exposed. 



B. The wave horn equation as a control system 
Suppose a dynamical system is given as : 



X ^fiX,U) (5) 



where X{x) G K" is the state, U{x) G [C/i, U2] C K the control and a prime still denotes the derivative with respect to the 
independent variable x. The control U is purposedly restricted to be a scalar as this suffices for the needs of the present 
work (see [23], [8] for a more general exposition). The purpose of Optimal Control theory is to find a control U such that 
some functional (a scalar function of the state and the control) : 

J{X, U)^ r g{X, U)dx (6) 

Jxo 

is maximized or minimized, making it an optimization problem under dynamical constraints, named for this reason infinite- 
dimensional. This variational-type problem is solved the following way : using a vector Lagrange multiplier /i, with 
components fj,i, i — 1, . . . ,n, define the Hamiltonian function of the problem, a scalar function, as : 

H{X,U,^l)=g{X,U) + ^l^f{X,U) (7) 

One then obtain necessary conditions for an optimum as the differential-algebraic system of equations (DAE) : 

X' ^ ^ 

- if (8) 



= 



dX 
dH 
dU 



together with suitable initial/boundary conditions. Notice first that, as is a scalar, all of the three above equations are 
vector ones : ^ e.g. is the vector with components i — 1, . . . ,n. One has then a set of 2n + 1 equations, the first 2n 
ones being differential equations, the last one being of order zero. Notice also that these conditions ([8]) are only necessary 
so that more investigations are needed to get sufficient conditions, i.e. a complete and unique solution. Also, the so-called 
Pontryagin maximum principle [23], [8] allows to have more precise results for the above when the control has e.g. bound 
constraints, as will be the case for the acoustic design problem here. 

Now, in order to fit the acoustic design problem at hand to this control theoretic setting, it is useful to reformulate the 
second order equation (O as a two dimensional first order system, suitable for control purposes in view. To this end, define 
the following variables : Xi — D,Y — <j), Z — (p {— v) and W the vector with components Xi, Y, Z. Also, as the section 
diameter D of the horn is an unknown to be determined as a function of x, it can be considered as a control variable to be 
designed. Actually, one can control either D itself or the way it varies along the x axis, that is one can control its derivative 
D with respect to x : as the modelling hypothesis i) above is that the section is a slowly varying function of x, one natural 
choice for the control is the derivative D and bounds will have to be imposed to it in order to satisfy the hypothesis given 
by inequality ([1]). As a consequence, defining the control variable U as D , constraints on it will be bound constraints : 
Di < D < D2- Together with those on the control, constraints can also be imposed to the state through design/building 
constraints on the diameter D itself : an obvious mandatory state constraint is e.g. Xi = D > 0. After having developped 
the second order derivatives, and with these notations, equation (O rewrites as the first order differential system : 

X[ = U 

y' = Z (9) 

Z' = -2^Z-PY 

Xi 

which is a dynamical control system, affine in the control U : W — f{W,U) = hi{W) + h2{W)U, with immediate 
definitions for hi and ft,2- One has even a linear drift : hi{W) — AW, where matrix A is clearly defined. The singularity in 
Xi = is not a real problem as it corresponds to a vanishing diameter, a highly uninteresting situation (except possibly at 
one boundary, the apex of a complete cone e.g.). Now, referring to theorem I, item 2 of the appendix, for each eigenvalue 



An = kf^ the eigenfunction (pn satisfies the wave horn equation with the same boundary conditions. One can then rewrite 
for each eigenvalue A„, the corresponding first order differential system of equations : 



^2n — X2n+1 



with : = ^n,X2n+i = V'n- This set of equations, i = 1,..., together with the equation ~ U constitutes the 
controlled dynamical model that will enter the optimization process as constraints. Remark that on a practical side, only a 
finite number of eigenvalues will be retained, as is detailed below in section |lll] One should notice also that, whereas in 
most control theoretic situations one looks for a feedback control, the searched after control for the above horn equation 
is intrinsically open-loop, as the duct is designed once for all (at least in the present state of technology...) so that optimal 
control theory and the Pontryagin Maximum Principle [23], [8] are well suited in the present context, although giving only 
necessary conditions, as already mentionned. 

III. Oscillation regime as an optimization criterion 

Besides the above simple, albeit new, problem rewriting, the crux in the followed approach of wind instrument design 
is in formulating the oscillation regime, as defined in section Ill-AI It must be done in a way suitable for conducting the 
optimization while having the controlled dynamical model of section III-BI above at hand. Considering the above definition 
for an oscillation regime, the potential (f> can be assumed to be a time-periodic function. For harmonicity requirements of 
the signals inside the duct, the wave numbers fc„s' are fixed data given as multiples of a fundamental frequency : fc„ = nko, 
ko — —fo- To each of the A:„s' corresponds a solution (pn{x) of the Webster equation with fc„ as wave number. The set 
{(Pn{x)}nei, is a basis of the Hilbert space of square-integrable functions on [0, L], thanks to Theorem 1, item 3, appendix. 
Thus, thanks to Theorem 1, item 4, appendix, the overall potential (p can be written as : 

which is seen to be the Fourier series decomposition of (l){x, t) with respect to the time variable. Thus Parseval theorem 
allows to write : 

The energy E can now be computed as a function of Cn,(pn,'Pn' simply substituting for (f> from ( ITTI ) into ([4|, which gives, 
thanks to the relations (fTSl i : 



T £ ^ + ^l^-n) = T ("^^ ^ -li\X2n+l\' + kl\X2n\')^ dx (13) 



E : 

the last equality stemming from the definition of X2„,X2„+i at the end of section HIl Now assume that the energy is 
conserved inside the horn. At least, one can consider in a first approximation that the small part of energy which is radiated 
outside the duct or dissipated at the inside boundary is exactly compensated with the energy brought in by the excitation 
mechanism such as the reed of the instrument. Such dissipation and other phenomena should surely be considered in future 
work, possibly through the boundary conditions. The definition of an oscillation regime leads then to maximize E as defined 
in equation ( fT3T l. Thus, with the ujns' known and fixed as given data, the optimization problem above has to find D{x) 
and the c„i'. Observe nevertheless that obtaining optimal c„s' does not give any indication on some perceptive quality 
factor which is indeed important for a high quality design. Thus it is likely at this stage that optimizing the above-defined 
oscillation regime is unsufficient to that end. Section [V] will present simulations for the above defined problem, after the 
design problem itself has been mathematically posed in the following section. 

IV. Optimal design of a horn shape 

Thanks to the results of section and |III1 the design of axisymmetric wind instruments can be now formulated as an 
optimal control problem in the following way. Consider an oscillation regime as defined and developped in sections Hi- Al HID 
each component being governed by the Webster horn equation with its own wave number fc„ and keeping the approximation 
of the solution to a finite number, N, of its first terms. The data are made of the set of fixed multiples of a fundamental 
frequency /o, leading to the set /c„ ~ jnko,jn — - ■■ ,N, with j„ allowing to take into account that a complete series 
of harmonics of fo or only the odd part of this series can be present. As cj) is unknown, the c„s' of its time-harmonic 
decomposition are unknown parameters and considered as constant controls. In the following, the c„s' will be gathered in 
a parameter vector noted C — (ci, C2, . . . , cat). Then the searched after shape of the duct and the c„s' are obtained as the 
solution of the following problem : 



(V) Maximize E as given by equation ( fT3] ). with respect to D and C 




subject to : 

1) The 2N + 1-dimensional dynamical model (n = 1, . . . , N) : 

= U 

= X2n+1 (14) 

each n*^ pair of the last 2N equations being a Webster equation with wave number A;„. 

2) Bound constraints on the variable control : Di < U = D < D2. 

3) State constraint : Xi > a > 0, a a fixed real number 

It is apparent that, whenever U > 0,Xi{0) > 0, the state constraint is satisfied. Within this context, the strong Pontryagin 
Maximum Principle (PMP) [23], [8] is applicable. Necessary conditions for an optimal control (i.e. an optimal duct shape) 
are obtained in the following way. Along the general setting given by (|5]l, (|6]l, O, (|8j, adjoin the dynamical constraints (fl4] i 
to the criterion E, (fT3T l. through Lagrange multipliers fi — {nn)n=i....,N , define the Hamiltonian of problem (V) as : 

H{X, U, C, fl) = ^Xf + klXi„) + filU + E1i(M2„^2„+1 - /i2„+l(2^X2„+i + klX2n)) ^^^^ 

Then, following (|8]i, jj. is the solution of the adjoint differential system (n = 1, . . . ,N) : 

— dH _ rrpo y v^Af „2 ^ -1^2 , 1.2 y2 \ n U ^ 

~ "axT ~ ^ri=l "^n(^2n+l + '«ri-^2nJ ~ ^IfT 2^„=l M2n+lA2„+l 

^^'2n = -^=kU^^2n+l~^clX^^X2n) (16) 

M2n+1 — ^ dX2„ + i ~ ^X;"M2n+l ^ fJ-2n ^C^A^ A2,i+1 

Thus the optimization problem will be complete when initial-boundary conditions are specified. Notice that time-initial 
conditions are not relevant here as the time variable has been eliminated through the hypothesis of time-harmonic regime 
and of energy conservation : the outputs of the optimization process are space-dependent, giving e.g. the mode shapes inside 
the duct. As a starting point and in order to simply illustrate the method, without going into deep physical considerations of 
a specific instrument, that are postponed to a forthcoming paper, the end at a; = is assumed to be closed (typically where 
the excitation would be placed) so that the velocity and thus (p vanish there : Vn = 1, . . . , N, (p2n{0) = -'^2n+i(0) = 0. 
Also the other end is assumed to be open (e.g. the other end of a brass instrument) so that, at the first order approximation, 
the pressure and thus (p vanish : Vn = 1, . . . , N, Lp2n{L) — X2n{L) = 0. These conditions lead to a Two-Point Boundary 
Value Problem (TPBVP) for the differential system (fT4l i.(fT6b : one part of the state is fixed at one end and another part is 
fixed at the other end. When some state component Xj is left unspecified at one end, the corresponding costate component 
/ij must vanish there [23]. Thus in the present situation : /i2n(0) = fi2n+i{L) = 0. Other conditions can be imposed with 
more realistic considerations. For computing the optimal control, the Pontryagin maximum principle implies that, for X 
(resp. fi) a solution of the state (resp. adjoint state) equation, an optimal control U and optimal parameter vector C are such 
that : 

H = H(X,U,C,fi) ^ max H(X,U,C,n) (17) 

{Di<U<D2],C=Cst} 

But one can observe that H is affine with respect to U so that : 

OH 2 ^ 

^ = Ml - X! ^2n+1^2n+l = (18) 

which is also named the switching function for this problem as its zeroes and its derivative can help determine the shape of 
the duct. As equation ( ITSl l does not allow to compute the control U explicitely, one is faced with a problem with singular 
extremal arcs [5], chap. 8. The standard method is to compute successive derivatives of dH/dU with respect to x until one 
is able to get an expression for U. In the present case, two such derivatives allow this, after tedious but straightforward 
computations. In case U is bounded, as it is the case here, the minimum of the hamiltonian H with respect to U, where 
dH/dU = 0, can occur at the boundary of the domain (see [5], [23]). This fact can be observed in the simulations : the 
control can be on its bound along several intervals of the integration interval (figure|4]l or can go from one bound to the other 
(figure |9]l, revealing a bang-bang type control [5]. The conclusion is that one generally obtains quasi-cones. Cones could not 
be obtained for exactly harmonic frequencies [7], chap. 7. This is compatible with the results obtained in [9]. Now, thanks 
to the developments of this section, one has sufficient material to illustrate the design method through a few simulations, on 
which some qualitative observations can be made and confirmation of the theoretical results is given. Remember nevertheless 



that in the following, only candidates for optimal shapes are obtained at this stage (necessary conditions) because, due to 
theorem 2, only q = 2D j D can be uniquely determined from the given boundary spectral data. Further data and theoretical 
investigations are necessary towards a really optimal shape. 

V. Numerical results 

The numerical results presented here must be considered as qualitative illustrations of the above theoretical results. No 
interpretation in terms of musical quality will be attempted at this stage. Such results will be pursued elsewhere more deeply 
in order to derive conclusions on realistic physical and perceptual basis. Nevertheless, some features are worth noticing at 
this stage, when faced with some characteristics of real instruments. 

A. Data 

In the following, data are as follows : the air mass density is taken to be po = 1, the sound velocity is that in free 
space at standard temperature : c — 3A0m/s and one focuses on a fundamental frequency for the note usually labelled A4, 
i.e. / — AAOHz, leading to the wavelength X = c/f — 0.772to. The chosen wave numbers are then indicated under each 
figure : fci — 27r/A, ki — riiki, i — 1, . . . ,n and Ui an even or odd integer. At the narrow end, x = 0, one has D{0) fixed 
with different values for each simulation (see below) and cf) (0) = 0. At a; = L, the condition 0(0) = is satisfied only 
approximately, by adjunction to the criterion through a simple penalty method. The duct diameter derivative is allowed to 
vary in the bounded interval [D^,Dj^j], D^,Dj^.j being indicated on each figure. They are chosen to satisfy inequality ([T]i 
with R — D/2, for the maximum value k of the fc„s', i.e. are within the validity domain of the wave horn equation, 
according to The interval of integration, i.e. the length of the duct, is taken to be L = A but could be taken as an 
unknown as well, leading to an analogous of the so-called time -optimal control. Last, all unspecified initial values for the 
state variables and all initial controls are taken to be random variables uniformly distributed on [0, 1]. 

B. Numerical method 

The numerical method used to solve the optimization problem exposed in section |IV] is the following (see e.g. [5], 
chap. 7 for more details on several possible numerical methods used for optimal control problems) : an initial guess is 
given for the control vector and, using this control, the state equation is numerically integrated from X{Q), considering 
too the unspecified initial conditions as constant controls, which are to be determined with the optimization process. The 
integration was performed thanks to a predictor-corrector scheme with an explicit Euler method for the prediction step and 
a Crank-Nicholson scheme for the correction. Then the costate equation is integrated backwards -notice that this equation 
is always linear in the costate- using the obtained terminal values for X{L). This allows to compute the gradient of the 
objective function with respect to the control, thanks to dH/dU [8]. Then a standard optimization routine is used, passing 
through the previous steps in a recursive way, in order to make this gradient decrease, until the specified terminal conditions 
at X = L are -approximately here- achieved and the criterion evolves x no more, to a specified precision (10^^). For the 
time being, the criterion and the terminal constraints have been gathered in a single objective function through a simple 
penalty method but this could be improved. For the optimization, a quasi-Newton method with projection has been used 
together with a BFGS method for an estimation of the hessian. Such an approach is known to have, as usually Newton-type 
methods, the drawback that it can give local minima and to be sensitive to the initial guess for the solution [10]. Thus, 
other methods such as direct solution of the TPBVP by shooting or multiple shooting techniques should be interesting to 
investigate but this is deferred to future work because an extra work has to be done in order to have an explicit expression 
for the control, as a function of the state and costate (see the discussion in section |IV] above). Nevertheless several different 
initial conditions for the state X were tested : the results, not exposed here, did not show to be that sensitive but this should 
be confirmed theoretically. 

C. Simulations 

The few simulations presented here have been done while fixing, respectively, two, five and ten components for the overall 
potential inside the duct, i.e. = 2, 5, 10 in equation (fT4] i. In each case, the duct shape is shown first, followed by the 
corresponding modes shapes that have been normalized to their maximum value, at x = 0. In addition for the case of two 
components only, the diameter derivative has been shown to illustrate the fact that it can be only piecewise continuous and 
that the phenomenon, mentionned in section HVl of singular extremal arcs joining regular arcs can appear (see figure |4|i : the 
control U = D is on its bound on some subintervals. For the five components cases, the shape is also quasi-conic but in a 
less obvious manner But for the ten component case, one has an example of abrupt change in conicity and D in figure |9] 
illustrates the possibility of bang -bang type control mentionned in section |IV] On each figure the imposed lower (Di) and 
upper (I?2) bounds are indicated. Notice that the impedance, p/u and the instantaneous power, pu {p, the pressure, u the 
volume velocity), for all x are easily obtained from the outputs, (p{x),(p {x) and D{x) Vx G [0,L], of the optimization 
process. The simulations are as follows : 

1) Two components : For this first simulation, the fundamental frequency is A4 and the second component is the second 
harmonic (double frequency). The diameter at the closed end is D{0) — 2cm. The duct shape is shown in figure |2l 



the two modal shapes in figure [3] D is shown in figure |4] as it makes appear the phenomenon described in section HVl 
where the control is at the bound on two subintervals, with the consequence that one obtains a quasi-cone for the 
duct, made of conic pieces joined by more complex but smooth parts. 

2) Five components : In this second simulation (see figures I5l6l l, the data are the same as for two components, except 
for D{0) = 1cm. The five components have frequencies fi = (i + 1) * /o; i = 0, . . . , 4. One has not shown D here 
but the result is here again a quasi-cone, although this is not apparent again on figure |5j because D varies much less 
but still varies along the duct. 

3) Ten components : In this third simulation (see figures I7I8I9I ). the data is the same as for five components, except for 
D{0) = 5cm. The ten components have frequencies fi — {i + 1) * fo;i = 0, . . . ,9. The results are noticeable as here 
again, three conic pieces are found so that the derivative D is shown in figure |9] making it here again a quasi-cone, 
with a bang-bang type behaviour for D . The duct is made of conic pieces with different conicity joined together. 
Only the five first modal shapes have been displayed, for better readability. 

D. Qualitative observations and commentaries 

The above simulations show only qualitative albeit important results at the present stage. The main one is that, as foreseen 
by the theoretical investigations, and for the specific chosen criterion, the duct shapes are piecewise continuous, this being 
particularly clear on the case of two and ten components ; observe too in both cases the fact that the control is at the 
bound on some subintervals, meaning that regular arcs and singular arcs can coexist. One can see also that quasi-cones are 
obtained, in an obvious way for the first and second simulations. Recall nevertheless that, as mentionned, the numerical 
method can be sensitive to initial conditions in the optimization process, so that it is difficult to interpret these results in a 
precise manner at the present stage. The qualitative observation is that the shapes are quasi-cones. At the same time, the 
modal shapes (i.e. the eigenfunctions of the Sturm-Liouville problem) show a behaviour that is qualitatively in agreement 
with what is theoretically predicted for cones (see [7], chap. 7), i.e. one has modal shapes that decrease as x increases to 
L. On one hand, these results are coherent with what is known for such geometries [9]. On another hand, the fact that one 
gets quasi-cones essentially would lend to conclude that the used notion of oscillation regime does not allow to capture 
all the important features that are sought after for musical quality, as flaring horns in brasses e.g. cannot be achieved with 
it : significant and important nonlinear phenomena [7] are not taken into account with the here retained model and design 
criterion. Thus, using more realistic models and refined criteria should be used instead. 

The computed modes have an extremum at the closed end, coherent with the chosen boundary condition (0 (0) = 0). At 
the other open end, the modes does not vanish exactly : one reason for this is linked to the way this end condition has been 
taken into account in the numerical method : a simple penalty method makes a compromise between satisfying the criterion 
and the boundary conditions there. The above suggested method of directly solving the TPBVP through shooting techniques 
should give better results from this point of view. 

VI. Conclusion and perspectives 

Thanks to a reformulation within a control theory framework, the design of axisymmetric wind instruments has been 
revisited. It has been shown that in order for an oscillation regime in the sense of A.H. Benade to take place inside the air 
column and with the considered linear model, the shape tends to be conic, as the number of time-harmonic components 
grows. This tends to show that such a criterion, as it has been mathematically described, is unsufficient to grasp the 
necessary conditions that lead to high-quality instruments. No other conclusion from this point of view is given at this stage. 
Nevertheless, and this was one main purpose of this work, the approach makes it very flexible to deal with a great variety of 
design constraints, either on the control or on the shape itself. At the present time, only qualitative results have been given 
to illustrate the theoretical results and the field is wide open to numerical investigations as well as to experiments on more 
physical and musical premisses. In that respect, some important issues can be straightforwardly put in perspective : 1) A first 
possible development is to use the method in a model fitting way, using an experimentally measured input impedance from a 
real instrument and taking a distance measure with the corresponding model output as criterion, in a similar way as in [18], 
[22]. More realistic duct shapes are to be expected and it will be interesting to compare results from different approaches. 
2) Refine the question of physical/perceptual criteria and of the horn model : this relies upon physical considerations as 
well as on experiments with real existing instruments played by expert musicians. This is surely a central point to get a 
practical tool useful for real design. One consequence of modifying the criterion is that the structure of the obtained shape 
will likely be different from that obtained in the present study. 3) Develop the method for woodwinds, thereby including 
toneholes in the design process. 4) Take into account some imposed-shaped parts such as pieces of cylinders or cones from 
built-in components, such as for brasses. This implies using path constraints in the optimal control framework. 5) All the 
previous points depend on the modelling question of how to include the nonlinear excitation mechanism (reed, lips...) in the 
model as well as other nonlinear phenomena, which were discarded in the present study. 6) Last, an interesting question is 
to try to retrieve, with the presented approach, the lattice of sound tubes results found in [9]. This will imply to modify the 
mathematical framework of the controlled dynamical system, i.e. the function spaces where the state and the control live, 



in order to be able to take into account discontinuous, and not only piecewise continuous, duct shapes. The mathematical 
theory of distributions could be an adequate framework for this. These items are but a small part of what should be addressed 
within the proposed framework to achieve a satisfying design methodology. 

VII. Acknowledgements 

The author is undebted to two anonymous reviewers an to the editor of the special issue, J. Kergormard, for their numerous 
advices that helped greatly improve the paper. 

References 



[1] Musical acoustics network. Summer meeting : wind instruments acoustics, h ttp://www.music.ed.ac.uk/euchmi/man/mxhta.litmll 2009. Edinburgh. 
[2] A.H. Benade. Fundamentals of musical acoustics. Dover, 1990. 2d edition. 

[3] H. Bouasse. Instruments a vent, I, II. Delagrave, 1929. New edition Librairie A. Blanchard, 1986. 

[4] M. Braneau. Manuel d'acoustique fondamentale. Etudes en mecanique des mateiiaux et des structures. Hermes, 1998. 

[5] A.E. Bryson and Y.C. Ho. Applied optimal control. Hemisphere Pub. Corp., 1975. Revised printing. 

[6] D.M. Campbell and C. Created. The musician's guide to acoustics. Oxford University Press, 1987. 

[7] A. Chaigne and J. Kergomard. Acoustique des instruments de musique. BeHn, 2008. 

[8] F.C. Clarke. Optimization and nonsmooth analysis. Wiley, New York, 1983. 

[9] J. P. Dalmont and J. Kergomard. Lattices of sound tubes with harmonically related eigenfrequencies. Acta Acustica, 2(5):421^30, 1994. 
[10] J.E. Dennis and R.B. Schnabel. Numerical methods for unconstrained optimization and nonlinear equations. Prentice Hall, Englewood Cliffs, New 
Jersey, 1983. 

[11] J. Dieudonne. Foundations of modern analysis. Academic Press, 1960. 

[12] E. Eisner. The complete solution of the "webster" horn equation. / Acoust. Soc. America, 41:1126-1146, 1967. 
[13] N.H. Fletcher and T. Rossing. The physics of musical instruments. Springer, 1991. 

[14] B.J. Forbes, E. Roy Pike, D.B. Sharp, and T. Aktosun. Inverse potential scattering in duct acoustics. J. Acoust. Soc. Am., 119(1):65— 73, 2006. 
[15] L. Henrique and J. Antunes. Optimal design and physical modelling of mallet percussion instruments. Acta Acustica u. with Acustica, 89:948-962, 
2003. 

[16] O. Inacio, L. Henrique, and J. Antunes. Design of duct cross sectional areas in brass-trapping resonators for control rooms. Noise Control Eng. J., 
55(2):172-182, 2007. 

[17] A. Katchalov, Y. Kurylev, and M. Lassas. Inverse boundary .spectral problems. Chapman & Hall/CRC, 2001. 

[18] W. Kausel. Optimization of brasswind instruments and its application in bore reconstruction. J. New Music Research, 30(l):69-82, 2001. 
[19] W. Kausel, P. Anglmayer, and G. Widholm. A computer program for optimization of brass instruments. In Forum Acusticum, Berlin, 1999. 
[20] G. Le Vey. Optimal control as a tool for the design of wind instruments. In Proc. International symposium on musical acoustics, ISMA07, Barcelona, 
sept. 2007. 

[21] G. Le Vey. Optimal control theory and design of wind instruments. In Proc. IEEE Mediterranean conference on control, MED08, Ajaccio, June 
2008. invited conference. 

[22] D. Noreland. Numerical Techniques for Acoustic Modelling and Design of Brass Wind Instruments. Acta universitatis upsaliensis, Uppsala University, 

2003. Comprehensive Summaries of Uppsala Dissertations from the Faculty of Science and Technology, 862. 
[23] L. Pontriaguine, V. Boltianski, R. Gamkrelidze, and E. Mitchenko. Theorie mathematique des processus optimaux. Mir, Moscow, 1978. 
[24] G.R. Putland. Modelling of horns and enclosures of loudspeakers. PhD thesis, Univ. Queensland, 1994. 
[25] S.W. Rienstra. Webster's horn equation revisited. SI AM J. Appl. Math., 65(6):1981-2004, 2005. 

Appendix 

a) The horn wave equation: Instead of the usual, Newton-type modelling approach, making use of forces and moments 
balance equations, a lagrangian one is followed here, as it fits better the variational context for the design methods that 
form the core of the present work. For standard derivation the reader is referred to [24], [4] for physical considerations and 
to [20] for control theory formulation. 

For a general pressure field, the lagrangian action density at each time instant inside the duct is the difference between a 
kinetic term and a potential term, that are computed in a standard way, as follows. Assumption ii) implies that the duct can 
be considered as a continuous stack of cross-sections S{x), parametrized by the abcissa x. For each section S{x) located at 
X along the horn axis, an action density is computed as the integral of densities of the particles over the section. This leads 
to an expression proportional to the cross-section area, i.e. to D'^(x). Firstly, the kinetic term writes : 



'Six) 

Similarly, the potential energy term is given as 
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Thus the lagrangian action density writes : 



L{x,t) = T{x,t) - U{x,t) = ^^m' - l-n (21) 
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Eventually, the lagrangian action inside the duct writes : 



C= [ L{x,t)dxdt^^ [ D^{\q^,\^ -\^\^)dxdt (22) 

Jo ° Jto Jo C 



According to Hamilton's stationary action principle, the dynamics inside the duct is obtained as the Euler equation of the 
above action £ : 

oL a , oL . o , oL , „ 

— — -A — ) = (23) 



But one can see that = 0, thus 



and eventually, dividing by 2D^ : 
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which is recognized to be the wave horn equation or Webster horn equation. The interest of the above variational derivation 
for this equation hes first in the fact that the optimal control approach to the design is variational in nature too. Thus, due 
to the harmonic nature of waves inside musical instruments, the horn equation writes, for one fixed value lu : 

2^'^=. + ^'Z' = (26) 
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i.e. using primes from now on to denote the derivatives with respect to the spatial dimension x along the axis, the only 
remaining independent variable : 

0"+2^(/.'+A:2,/, = O (27) 
this equation being adjuncted a set of suitable boundary conditions : 

ai(/)(0) + (0) = , a2<j){L) + fea^' [L) = (28) 

which represents a resonator without active components and for which losses can be taken into account in the boundary 
conditions that will be precised for the simulations in section [V] Thus one is faced with a homogeneous Sturm-Liouville 
problem. 

b) Elements of Sturm-Liouville theory: Two theorems from spectral theory of Sturm-Liouville problems are recalled 
here for self-containedness : one [11] for the so-called "direct problems" and the second [17] for "inverse problems". Let 

/ = [O, L] and q{x) = 2^. Then : 

Theorem 1 [11] : For every function q{x) continuous in / : 

1) The Sturm-Liouville problem has an infinite strictly increasing sequence of eigenvalues A„ G M such that lim„_j.oo A„ = 
+CXD and the series converges. 

2) For each eigenvalue A„, the homogeneous Sturm-Liouville problem has a real-valued solution tfn{x) such that 

(p'^{x)dx — 1, which is unique up to a multiplicative real constant. 

3) The sequence {(fin) is an orthonormal system in a convenient Hilbert space of functions. 

4) Let w be complex-valued continuous function defined in /, the primitive of a ruled function w such that : (i) w 
is continuous in /, except possibly at a finite number of interior points, (ii) w has a derivative w continuous in 



every interval where w is continuous, (iii) w satisfies the boundary conditions in ( l28l l. Then, if c„ w, ipn >= 
w{s)(pn{s)ds, one has : w{x) = J^n'^n'Pn where the series converges uniformly and absolutely in /. 
On another hand, the design question itself relies upon the following inverse boundary spectral theorem : 
Theorem 2 [17] : Assume that {Ai, A2, . . . , ^fiiiO), ^2(0), ■ ■ ■} are the boundary spectral data of the Dirichlet-Schrodinger 
operator, ^0 = + 9' corresponding to the above Sturm-Liouville problem, on an interval [0,L]. Then, these data 

determine L and q{x) uniquely. 

It is moreover worth noting that only q can be determined uniquely from these data. This means here that only the ratio 
D /D is so. Thus to determine D itself needs supplementary data. 
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Fig. 3. Modal shapes, 2 components, —0.2 < D < 0.2 
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Fig. 5. Duct shape, 5 components, —0.2 < D < 0.2 
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Fig. 9. Derivative of the duct diameter, 10 components, —0.2 < D < 0.2 



